#!/usr/bin/env python
#-*- coding:utf-8 -*-

from math import *
from sys import *


# como a matriz gerada eh ortonormal, a inversa dela eh a transposta!
# simples.



def vetor_norma(vet):

	n = len(vet)

	a = 0
	for i in range(n):
		a += vet[i] * vet[i]

	return sqrt( a )


def vetor_prod_escalar(v1, v2):

	n, n2 = len(v1), len(v2)

	if n != n2:
		raise Exception( "Os vetores não possuem o mesmo tamanho!" )

	a = 0
	for i in range(n):
		a += v1[i] * v2[i]

	return a


def vetor_mul(vet, s):

	nvet = []
	for i in range( len(vet) ):
		nvet.append( vet[i] * s )

	return nvet


def vetor_add(v1, v2):

	n, n2 = len(v1), len(v2)

	if n != n2:
		raise Exception( "Os vetores não possuem o mesmo tamanho!" )

	nvet = []
	for i in range(n):
		nvet.append( v1[i] + v2[i] )

	return nvet


def print_vetor(nome, vet):

	print( "Vetor " + nome + ":" )

	n = len(vet)

	for i in range(n):
		print( "   " + str( vet[i] ) )

	print( "   Norma: " + str( vetor_norma(vet) ) )
	print( "" )



def print_matriz(nome, mat):

	print( "Matriz " + nome + ":" )

	n = len(mat)

	for i in range(n):
		for j in range(n):
			stdout.write( "   " + str(mat[i][j]) + " " )
		stdout.write("\n")

	print( "   Norma das linhas:" )
	for i in range(n):
		print( "   Linha " + str(i) + ": " + str( vetor_norma( mat[i] ) ) )

	print( "   Produto Escalar entre as linhas:" )
	for i in range(n):
		for j in range(n):
			if i != j:
				print( "   Linha " + str(i) + " com " + str(j) + ": " + \
					str( vetor_prod_escalar( mat[i], mat[j] ) ) \
				)
	print( "" )



def matriz_x_vetor(mat, vet):

	n = len(vet)
	resp = [0.0] * n

	for i in range(n):
		for j in range(n):
			resp[i] += mat[i][j] * vet[j]

	return resp


def matriz_transposta(mat):

	n = len(mat)
	tmat = matriz_n_x_n( n )

	for i in range(n):
		for j in range(n):
			tmat[i][j] = mat[j][i]

	return tmat


def matriz_n_x_n(n):

	mat = []
	for i in range(n):
		mat.append( [0] * n )

	return mat


def dct(x):

	n = len(x)
	mat = matriz_n_x_n( n )

	for i in range(n):
		mat[0][i] = 1.0 / sqrt(n)

	for i in range(1, n): # de 1 ate n
		for j in range(n):
			mat[i][j] = \
				sqrt( 2.0 / n) * \
				cos( i * pi/(2*n) * (2*j + 1) )

	print_matriz( "DCT", mat )
	return mat



def idct(x):

	n = len(x)
	mat = matriz_n_x_n( n )

	for i in range(n):
		mat[0][i] = 1.0 / sqrt( n ) # no original, havia isso, que eh igual: sqrt( 1.0 / n )

	for i in range(1, n):
		for j in range(n):
			mat[i][j] = \
				sqrt( 2.0 / n ) * \
				cos( i * pi/(2*n) * (2*j + 1) )

	print_matriz( "IDCT", mat )

	tmat = matriz_transposta( mat )
	print_matriz( "IDCT transposta", tmat )

	return tmat




def main():
	entrada = [1, 2, 3]

	print_vetor( "Entrada", entrada )
	m = dct( entrada )
	saida = matriz_x_vetor( m, entrada )
	print_vetor( "Saida", saida )

	somatorio = saida[0] / m[0][0]
	print( "Somatorio do vetor de entrada: " + str(somatorio) )

	media = somatorio / len(entrada)
	print( "Media do vetor de entrada: " + str(media) )

	df = saida[2] / m[2][0]
	print( "Derivada segunda por diferencas finitas (?): " + str(df) )

	# na matriz 3x3, a linha0 contem: 0.57, 0.57, 0.57.
	# entao, ao fazer saida[0] + n, sera adicionado n x 0.57
	# em todos os componentes igualmente.

	#saida[0] += 1.0

	# se quisermos somar X, entao: X = n x 0.57 => n = X / 0.57.
	# assim, pra adicionar "10" em todo mundo:

	#saida[0] += 10.0 / m[0][0]

	# pra somar 10 em todos os componentes
	# saida[0] += 10.0 / m[0][0]

	# pra multiplicar por 2 todos os componentes
	#saida = vetor_mul( saida, 2 )

	# para somar a [media dos componentes originais] em todos os componentes
	#saida[0] = saida[0] * 2

	# para somar a [2 x media dos componentes originais] em todos os componentes
	#saida[0] = saida[0] * (1.0 + 2)

	# para subtrair a [3 x media dos componentes originais] em todos os componentes
	#saida[0] = saida[0] * (1.0 - 3)

	# resumo:
	# adicionar em saida[0] -> adicionar fator igual em todos
	# multiplicar em saida[0] -> adicionar n.media em todos

	# elemento da matriz [1][0] = 0.707106781187 = cos(45)
	# linha
	#saida[2] = saida[2] + 1.0


	print_vetor( "Saida modificada", saida )
	m = idct( saida )
	recuperado = matriz_x_vetor( m, saida )
	print_vetor( "Recuperado", recuperado )


if __name__ == '__main__':
	main()
